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Abstract 

Recently it has been observed that the bivariate generalized linear failure rate distribu¬ 
tion can be used quite effectively to analyze lifetime data in two dimensions. This paper 
introduces a more general class of bivariate distributions. We refer to this new class of 
distributions as bivariate generalized linear failure rate power series model. This new class 
of bivariate distributions contains several lifetime models such as: generalized linear failure 
rate-power series, bivariate generalized linear failure rate and bivariate generalized linear 
failure rate geometric distributions as special cases among others. The construction and 
characteristics of the proposed bivariate distribution are presented along with estimation 
procedures for the model parameters based on maximum likelihood. The marginal and 
conditional laws are also studied. We present an application to the real data set where our 
model provides a better fit than other models. 

Keywords: Bivariate generalized linear failure rate distribution; EM algorithm; Maximum like¬ 
lihood estimator; Power series class of distributions. 
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1 Introduction 


Compounding continuous with discrete distributions have been introduced and studied in the 


recent years. These method allow us to obtain distributions with great flexi 
to dev elop more realistic statistical models in a great variety of applications. 


jility and are useful 


Marshall and Olkin 


(| 19971 ) introduced a class of distributions which can be obtained by minimum and maximum 


of independent and identically distributed c ontin uous random variables, where the sample size 
follows geometric distribution. 


Silva et al. 


(|2013b introduced a class of distributions obtained 
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by mixing extended Weibull and power series distributions and stu died several of its statistical 


properties. This class cont ains the Weibull-geometric distribution (jMarshall and Olkin . 


Barreto-Souza et al 


to 


Silva et al 


(|2oi^ l for a brief literature review about some univariate distributions obtained 


1997 


20111 ) and other lifetime models studied recently. The reader is referred 


by compounding. 

The t hree-parameter gen e ralize d linear failure rate (GLFR) distribution has been intro¬ 
duced bv Isarhan and Kundu ( 2009 1. The hazard function of GLFR distribution can be in¬ 
creasing, decreasing and bathtub shaped, and it has the following cumulative distribution 
function (cdf) and probability density function (pdf), respectively, 


FQ{x]a) = {l — e X > 0, q;,/3,7>0, (1.1) 

/G(x;a) = a(;S + 7x)e-^^-2^'(l-e-^^-2^')"-^ (1.2) 


Recently, many studies have been done on GLFR distribution, and so me authors have 


extended it: the generalized linear ex ponential (jMahmoud and Alam 


201 


201C), beta-linear fail- 


( Jamkhanehl. 

2014 

1, Mel 

Donald-GLFR (1 

llbatal et al. 

2niJl. Poisson-G 

LFR ( 

Gordeiro et al. 


2014 1. GLFR-geometric (jNadaraiah et al. 


20141 1 and GLFR-power series (jAlamatsaz and ShamsI . 


1 are some univaria te extension of GLFR distribution. 


Sarhan et al 


(1201 il l have introduced a new bivariate distribution using the GLFR distribu¬ 


tion and derived several interesting properties of this new bivariate distribution. The proposed 
bivariate GLFR (B GLFR) distribution is an extension of bivariate generalized exponential 
(BGE) distribution (jKundu and Guptal. l2009l ) and its cdf, the joint pdf and the joint survival 
distribution function are in closed forms. This bivariate distribution has both an absolute con¬ 
tinuous part and a singular part, and is extended to a multivariate distribution. The cdf of the 
BGLFR model is given by 


-^bg(3:i, X 2 ; ai, Q!2, as) = 


if Xi<X2 


if Xi > X2- 


(1.3) 


In this paper, we compound the BGLFR distribution and power series class of distribu¬ 
tions and define a new class of bivariate distributions. This class contains the BGLFR and 
GLFR-power series (GLFRPS) distributions and is called the bivariate GLFR-power series 
(BGLFRPS) class of distributions. This paper is organized as follows. In section [2l we intro¬ 
duce the BGLFRPS distributions and obtain some properties of this new family. Some special 












































































models are studied in detail in Section [3l An EM algorithm is proposed to estimate the model 
parameters in Section [H A real data application of the BGLFRPS distributions is illustrated 
in Section [H 

2 The BGLFRPS class 


Let the random variable be a discrete random variable having a power series distribution 
(truncated at zero) with probability mass function (pmf) 


P{N = n) = 


an0^ 

C(0)’ 


n = 1,2, 


( 2 . 1 ) 


where > 0 depends only on n, C(0) = 0 ^ (0;®) (s can be oo) is such 

that C{0) is finite. Table [T] summarizes some particular cases of the truncated (at zero) power 
series distributions (geometric, Poisson, logarithmic, bi nomial and neg ative binomial). Detailed 


NoackI (|1950I L Here, C'{0), C''{0) and 


properties of power series distribution can be found in 
C"'{0) denote the hrst, second and third derivatives of C{0) with respect to 0, respectively. 


Table 1: Useful quantities for some power series distributions. 


Distribution 

O-n 

C{0) 


C'{0) 

C"i0) 

C"'(0) 

S 

Geometric 

1 

0{l-0)- 

-1 

(1-0)-^ 

2{l-0)-^ 

6(1-0)-^ 

1 

Poisson 

n\~^ 

e^-1 



e® 


OO 

Logarithmic 

n~^ 

- log(l - 

0) 

(l-0)-i 

(l-0)-2 

2(1 - 0)-^ 

1 

Binomial 

0 

(1 + 0)^- 

- 1 

k 

(e+i)!-* 

k{k-l) 

(e+i)^-* 

k{k-l){k-2) 

00 

Negative Binomial 

/n-l\ 

Qk 



k(k+ 2 e-l) 


1 

U-J 

(l-6»)'= 


(i-0)k+i 

6|2-fc(l_0)fc+2 

03-fc(l_0)fe+3 



Now, suppose {{Xin, X 2 n)',n = 1,2...} is a sequence of independent and identically 
distributed non-negative bivariate random variables with common joint cdf Tx(.,.), where 
X = (Xi,X 2 y. Take to be a power series random variable independent of {Xu, X 2 i). Let 

Tj = max{Aii,... AjAf}, i = l,2. 

For given N = n, the joint cdf of = (Y}, Y 2 ) is 

^Yi,Y2|Ar(yi,2/2|n) = iFx{yi,y2)T- (2.2) 

Therefore, the joint cdf of Y becomes 

TP ^ ^ ^fTp T C{0Fx{yi,y2)) 

FY{yi,y2) = 2_^{Fx{yi,y2)) =- QTf\ -• 

n=l '' ' '' ' 


(2.3) 


















In this case, we call Y has a bivariate F-power series (BFPS) distribution. 
The corresponding marginal cdf of Yi is 

c{eFxM)) 


FrAVi) = 


m 


i = l,2, 


and recently this univaria te class is considered 


expo nential-power serie s (jMahmoudi and Jafari 


ries (Flores et al 


201J), 


3 y ma ny authors: for example the generalized 


2012l b complementary e xponential-power se- 


complementary extended W eibull-power series (jCordeiro and Silval . 


2014) and and GLFRPS ( Alamatsaz and Shams . 2014) distributions. 


In this paper, we take F to be the bivariate GLFR distribution given in (jl.Op . Therefore, 
we consider the BGLFRPS class of distributions which is defined by the following cdf: 


FYiyi,y2) = < 


C ( “ 7^1 )“1+“3 - 5!/2 )“2 ) 

^- m - - if 2/1 ^ 2/2 

cm 


if yi > 1 / 2 , 


{ C{0Fa{yi;ai+a3)Fa{y2-,a2)) 
0(8) 

CidFa{yi-,ai)FG(y2-,a2+a3)) 

0 ( 8 ) 

and is denoted by BGLFRPS (ai, 02 , 03 ,/S, 7,6*)- 


if yi < 1/2 
if yi > 1 / 2 , 


(2.4) 


Theorem 2.1. Let Y = ( 11 , 1 ^ 2 ) has a BGLFRPS (oi, 02 , 03 ,/3, 7 ,0) distributions. Then the 
joint pdf of Y is 


/i(yi,y2) if 0 < yi < ?/2 

fY{yi,y2) = { f2{yi,y2) if 0 < y 2 < yi 

/o (y) if 0 < yi = y 2 = y, 


(2.5) 


where 


/i(yi,y 2 ) = 


/ 2 (yi,y 2 ) = 


/o(y) = 


C{0) 


/G(yi; ai + otf)fG{y 2 ] 02) [ 0 FG(yi; ai -|- aA)FQ{y2\ 02) 


( 2 . 6 ) 


xC"'(6'FG(yi;ai + oi'i)FG{y2]a2)) + C'{6FQ{yi;ai -h a3)-/^G(y2; 0 : 2 ))] , 

a 

/G(yi; ai)/G(y 2 ; a 2 + a 3 ) [ 6 '-/^G(yi; ai)^G(y 2 ; 02 + 03) 


C{0) 


(2.7) 


xC'"(6'FG(yi;ai)^G(y2;a2+ 0 : 3 )) +C"(6'-/^G(yi;ai)-/^G(y2;a2 + 0 : 3 ))], 

603 


C’{d){oii (y.2 T 03 ) 


foiy; ai + a2 + a3)C'{9FG{y; ai + a2 + as)). ( 2 . 8 ) 


Proof. It is obvious. 


□ 


As a special case, we cons ider C{9) = 6 + 6^^ (see also 


Morals and Barreto-Souzal . 


Mahmoudi and Jafaril . 


2012 


201l|). The pdf of the BGLFRPS class of distributions are depicted 


in Figure [U for /3 = 7 = 1 and some values of other parameters. 






































Figure 1: The pdf of the BGLFRPS class of distribution for some values of parameters: ai = 
a 2 = = 1, 9 = 1 (left top), ai = 02 = as = 1, 0 = 2 (right top), ai = 02 = as = 1, 0 = 0.5 

(left bottom), ai = 02 = 0.5, as = 1, 0 = 1 (right bottom). 


Proposition 1. Let (1^,12) has a BGLFRPS (ai, 02, as,/ 3 ,7, 0 ) distributions Then 

1. Each Yi has a GLFRPS distributions with parameters a^ + as, /3,7 and 9 with following cdf: 


F{x) 



C{9) 


X > 0. 


2. The random variable U = max(Y’i,l 2 ) has a GLFRPS distributions with parameters ai + 
0:2 + «Sj /5,7 and 9. 


3. If C{9) = 9, then Y has a BGLFR distribution with parameters ai, a 2 , as, /? and 7. 


4. p(yi < Y2) = 


Ql 


Oil+Q:2+Ct3 ' 


Proposition 2. Let iV(2/1 ? 2/2) cdf of BGLFRPS distributions given in (| 2 . 4 p . Then 

00 

Pr(?/i,y2) = ^PnPBG(2/i,2/2;nai,na2,nas), 

n=l 








where pn = P {N = n) = Therefore 

OO 

fi{yi,y 2 ) = '^PnfG{yi-,nai Pnas) fG{y 2 ;na 2 ), 

n=l 

OO 

f2{yi,y2) = '^PnfG{yr,nai) fG{y2]na2 + na3), 

n=l 

OO 

fo{y) = -^^-y^Pn/c (y;nai + na 2 + nas), 

ai + 02 + as , 

n=l 

where /g (•;««) is the pdf of GLFR distribution with parameters na, fd and 7. Note that 
/g (.; na) is the pdf of random variable max {Ui, ..., Un) where Ui ’s are independent random 
variables from a GLFR distribution with parameters a, /3 and 7. 


Proposition 3. The joint pdf of the BGLFRPS distributions provided in Theorem \2.1\ can be 
written as 

fY{yi,y2) = — ^ — 50(^1, ^2) H- 7 -^. — 


where 


9a{yi,y2) = 


e 


ai + 0-2 + 03 ai + 02 + 03 

01+02 + 03 / /i (yi, 2 / 2 ) if 2/1 < 2/2 


oi + 02 


/2 ( 2 / 1 , 2 / 2 ) if 2/2 < 2 / 1 , 


5 s( 2 /) = 777 ^/ 0 ( 2 /; ai+ 02 + 03 )C"( 6 IFg( 2 /;oi+ 02 + 03 )) if yi=y 2 =y, 

C[d) 

and 0 otherwise. Glearly gai-,-) is the absolute continuous part and gs{.) is the singular part. 
If as = 0 , it does not have any singular part and it becomes an absolute continuous pdf. 


Proposition 4. The conditional distribution of Yi given I 2 < 2/2 is an absolute continuous 
distribution with the following cdf: 


P(Yi < yi\Y 2 < 2 / 2 ) 




if 2/1 < 2/2 
if 2/2 < 2 / 1 - 


Proposition 5. The limiting distribution of BGLFRPS when 0 —>■ O'*" is 


lim FY{yi,y 2 ) 
e^o+ 


C{9Fx{yi,y2)) 

— m — 

aciFx{yi,y2)r + En=c+i(^nO^-%l^x{yi,y2)r 

0 ™+ «c + E“=c+l«n 0 "-'= 














{Fx{yi,y2)y 


{ (1 y^<y^ 

= I (i _,-/>«-?»;)“■ (i _ ,-te-i,iy<»+“>> if 

which is the pdf of a BGLFR distribution with parameters cai, ca2, caa, f5 and 7, where 
c = min {n G N : a„, > 0 }. 


Based on (j 2 . 2 jl . {Yi,Y2) given N = n has a BGLFR with parameters nai, na2, nas, (3 and 


7, and therefore, the joint pdf of (Yi, 12)-^) is 


fYi,Y2,N{yi,y2,n) 


' ^fin(yi,y2) 
< ^f2niyi,y2) 
. ^^fon{y) 


if yi < y2 
if y2 < yi 
if yi = y2 = y, 


where 


fin{yi,y2) 


/ 2 n( 2 / 1 , 2 / 2 ) 


fon{y) 


n^(,d + 7 yi )(,5 + 72 / 2 )(ai + 03)026 livi+y^) 

_ g-^!/l-i!/i^^(«l+«3)-l('2 _ Q-hy2-^y2')na2-l^ 
n ^(/3 + 72 /i)(/? + 72 / 2 )(a 2 + 03 )oie“^^^i+^ 2 ^“^ 

X(1 _ _ ^-Py2-^yl'^n{a2+a3)-l ^ 

n (/3 + -fy)a3e-^y-2y\l - e-/3y-iJ/")’^(«i+«2+“3)-i, 


The conditional pmf of N given Yj = yi and Y2 = 2/2 is 


/iV|Yi,Y 2 (?^l 2 /l, 2 / 2 ) 


n?an{0Ai{yi,y2))"~^ 

ki{yi,y2) 

rF an(6A2(yi ,y2)Y~^ 

k2{yi,y2) 

na„(dAo(y})"~^ 

' ko(y) 


where 


if 2/1 < 2/2 

if 2/2 < 2/1 

if 2/1 = 2/2 = y, 


( 2 . 9 ) 


ki{yi,y2) = 0FG{yi]ai + a3)FG{y2]a2)C"{9FG{yi-,ai + a3)FGiy2]a2)) 
+C"(6'Fg(2/i; «i + a3)FGiy2] 02)), 

^2(2/1,2/2) = 6'FG(2/i;ai)-?^G(2/2;a2 + a3)C'"(^'-^G(2/i;ai)-?^G(2/2;a2 + as)) 

+C'(6'Fg(2/i; ai)FG(2/2; a2 + 03)), 
ko{y) = C'(6'FG(2/;ai + 02 + 03)), 

(2/1,2/2) = (1 - - e-Py2-U)-2 = 


and 










^2 (2/i, 2/2) = (1 - e _ e lylY'^+^^ = FQ{yi]ai) FG{y2]a2 +a^) 

^0 [y) = (1 - Fg {y; ai + aa + as) 


Since e^C"'{e) + WC” {9) + C {9) = J2n=i n^an9^-\ 0C"{9) + C (9) = i n^an9^~^ 


and C (9) = nan9^ we can obtain the conditional expectation of N given Yi = yi and 


Y 2 = y 2 as 


E{N\Yi,Y2) 




Bi{yi,y2) 

ki{yi,y2) 

B2{yi,y2) 

k2{yi,y2) 

eAo{y)C"{eAoiy))+C'{eAoiy)) 

ko{y) 


if yi < y2 

if y2 > yi 

if 2/1 = 2/2 = y- 


( 2 . 10 ) 


where 


S,(yi,2/2) = {9My^,y2)fC"'{9Ai{y^,y2)) 
+WAi{yi,y2)C''{9Ai{yGy2)) 

+C'{9Ai{yi,y2)), z = l,2. 


Remark 2.1. If we consider Zi = min {Xji,... Xj^v}, i = 1,2, another class of bivariate 
distribution is obtained with the following joint cumulative survival function: 


FzuZ2izi,Z2) = P{Zi > Zi,Z2 > Z2) = 


C{9Fx{zuZ2)) 

C{9) 


where Fx{xi,X 2 ) = P{Xi > xi,X 2 > X 2 ). This uni variate class of distributions is studied 


in literature: for exa mple the exponential-power series 1 Chahkandi and Ganiali. 


Weib utll-yower series iMorais and Barreto-Souza 


201 A). Burr XII-power series IlSilva and Cordeird. 


200f hamltlw 


Bourauianon et al. 

2014 

) and linear failure rate-power series 1 

Mahmoudi and Jafan 


201 A) distributions. 


3 Special Cases 

In this section, we consider some special cases of BGLFRPS distributions. 

3.1 Bivariate GLFR-geometric distribution 

When C (9) = (0 < 9 < 1), the power series distribution becomes the geometric distribution 

(truncated at zero). Therefore, the cdf of bivariate GLFR-geometric (BGLFRG) distribution 




































is given by 


FY{yi,y2) = < 




if 2/1 < 2/2 
if 2/1 > 2 / 2 , 


and its pdf is given in (|2.5p with 

/i( 2 /i, 2 / 2 ) = (1 - 6 ')/G( 2 /i;ai + a3)/G(2/2;a2) 

72 ( 2 / 1 , 2 / 2 ) = (1 - 6 ')/G( 2 /i;ai)/G( 2 / 2 ;a 2 + a 3 ) 

fo{y) = 


1 + OFciyiiai + «3)-^g(2/2; « 2 ) 

(1 - 6'Fg(2/i; ai + a3)FG{y2] 02 ))^ ’ 
1 + 6>Fg(2/i; «i)-Fg(2/2; «2 + «3) 

(1 - 6'Fg(2/i; ai)-FG(2/2; 02 + a3))^ ’ 

(1 - 6 )a 3 fG{y; ai + 0:2 + 0:3) 


(ai +02 + q; 3)(1 - (^FG{y;ai +02 + 03))" 
Remark 3.1. When 9* = 1 — 9, we have 


( 2 / 1 , 2 / 2 ) = < 


0 * )“\l-e“'^^ 2 " 5 *'2 


if 2/1 < 2/2 
if 2/1 > 2 / 2 - 


/2 is also the cdf for all 9* > 0 fsee lMarshall and Olkin . 11991 ). In fact, this is in Marshal- 
Olkin bivariate class of di stributions. Also , the m arginal distribution of Yi is GLFR-geometric 
distribution introduced by 


Nadaraiah et al. 


mu)- 


3.2 Bivariate GLFR-Poisson distribution 

When On = ^ and C {9) = e^ — 1 (0 > 0), the power series distribution becomes the Poisson 
distribution (truncated at zero). Therefore, the cdf of bivariate GLFR-Poisson (BGLFRP) 
distribution is given by 


Fy ( 2 / 1 , 2 / 2 ) 
and its pdf is 




exp|e(l-e-'’''l-^^l)“\l-e-^«2-^!/i)“2+“3|_3 

e"-! 


if 2/1 < 2/2 
if 2/1 > 2 / 2 - 


/l ( 2 / 1 , 2 / 2 ) 


f2(91,92) 


fo{y) 


Ofciy+ai + 0:3) f0(92; 0:2) ew{OFG{yi-,ai + a3)FG(2/2; 02) - 1} 

X [ 9 FGiyi;ai + 03)^0(2/2; a 2 ) + 1 ], 

6fG(y+,ai)fG{y2-,a2 + 03) ex.p{ 9 FG{yr, ai)FG(y2; a2 + (13) - 1} 
X [6'Fg(2/i; ai)-FG(2/2; 02 + 03 ) + 1 ] , 

(q!1 -|- 02 + 03) 


/g (2/; 01 + 02 + 03) exp{ 6 »FG (2/; ai + 02 + 03) - !}• 





















3.3 Bivariate GLFR-binomial distribution 


When ttn = (^) and C (0) = (0 + 1)^ — 1 (0 > 0), where k {n < k) is the number of replicas, 
the power series distribution becomes the binomial distribution (truncated at zero). Therefore, 
the cdf of bivariate GLFR-binomial (BGLFRB) distribution is given by 


Fv (2/1,2/2) = < 






if 2/1 < 2/2 
if 2/1 > 2/2, 


and its pdf is 
/i (2/1,2/2) = 

/2 (2/1,2/2) = 

k{y) = 


k0 


/g(2/i; «i + a3)/G(2/2; a2)[0FG{yi;ai -|- q;3)Fg(2/2; 012 ) + 1] 


k-2 


(0 + 1)^-! 

X [/c6»Fg(2/i; «! + ak)FG{y2] 02 ) + 1] , 
k0 

- - -z - /g(2/i; ai)/G(2/2; 02 + a3)[0FG{yi; ai)TG(2/2; 02 + 03 ) + 1] 

(0 + 1 ) - 1 

X [kOFGiyi] ai)FG{y2; 02 + 03 ) + 1], 
k 0 a 3 fG { y;ai -ba2 + 03) 


k-2 


[(0 + 1) — l](oi -|- C 12 T 03 ) 

3.4 Bivariate GLFR-logarithmic distribution 


( 6 'Fg( 2 /; Oi + 02 + 03) + 1 ) 


fc-i 


When On = ^ and C{0) = —log (1 — 0) (0 < 0 < 1), the power series distribution becomes the 
logarithmic distribution (truncated at zero). Therefore, the cdf of bivariate GLFR-logarithmic 
(BGLFRL) distribution is given by 


Fy (2/1,2/2) 


and its pdf is 


log(^l-6»(l-e-'’^i-7S'i)“^’^“®(l-e-^5'2-5a2)“^j 

log(l-6») 

log^l-6»(l-e“^^l~T*'l)°’^(l-e“^*'2-T!/2)“^’'’“® j 
. log(l-6») 


if 2/1 < 2/2 
if 2/1 > 2/2, 


71 ( 2 / 1 , 2 / 2 ) 

72(2/1,2/2) 

fo{y) 


_ -^ 7 g( 2 /i; oi + 03)70(2/2; 02) _ 

log(l - 0)(1 - 0Fg(2/i; oi -b 03 )Fg(2/2; 02))^ ’ 

_ -^ 7 g( 2 /i; oi) 7 g( 2 / 2 ; 02 + 03) _ 

log(l - 0)(1 - 0 Fg(2/i;oi)Fg(2/2;o 2 -b 03))^’ 

_ - 0 O 37 g( 2 /; Ol + 02 + 03) _ 

log(l - 0 )(q;i -b 02 + 03) (1 - 0Fg(2/;oi -b 02 + 03))' 


3.5 Bivariate GLFR - negative binomial distribution 

When an = (^Zi) and C{0) = (yz^)^ (0 < 0 < 1), the power series distribution becomes 
the negative binomial distribution (truncated at zero). Therefore, the cdf of bivariate GLFR- 



if yi < y 2 


negative binomial (BGLFRNB) distribution is given by 


FY{yi,y2) 




if 2/1 > 2/2, 


and its pdf is 


/i(yi, 2 / 2 ) 


/2 ( 2 / 1 , 2 / 2 ) 


/o(2/) 


fc(l - 6 *)^/g(2/i;« 1 +«3)/G(2/2;«2)i^G~^(2/i;«l +«3 )^g~H2/2;«2) 

(1 - 6'Fg(2/i; ai + a3)FG(2/2; 02 ))^’^^ 

X [fe + 0 Fg(2/i; ai + a2,)FQ{y2] 02 )], 

fc(l - 6 *)^/g(2/i; «i)/g( 2/2; «2 + a3)-^G~^(2/i; «i)-^g~^(^2; «2 + «3) 
(1 - 6'Fg(2/i; a:i)FG(2/2; 02 + 03 ))^"^^ 

X [A: + 6'Fg(2/i; ai)FG(2/2; 02 + as)], 

fca3(l - effciyiai + a2 + a3)FQ~^{y;ai + a2 + as) 

(ai + 02 + a3)(l - OFciy; ai + 02 + 03))*"+^ 


4 Estimation 


In this section, we consider the estimation of the unknown parameters of the BGLFRPS dis¬ 
tributions. Let ( 2 / 11 , 2 / 12 ), • • •, ( 2 /mi, 2 /m 2 ) be an observed sample with size m from BGLFRPS 
distributions with parameters © = (ai, 02 , a 3 ,/3, 7 ,0)^ Also, consider 

Iq = {i ■ yii = 2/2j = yi} , h = {i ■ yu < 2/2i} , L2 = {z : z/ii > 2/2i} , 


and 

mo = |/o|, "ii = |Li|, m2 = 1/2!, m = mo + mi + m2 

Therefore, the log-likelihood function can be written as 

m = E log(/o(2/i)) + log(/l(2/li,2/2i)) + ^log(/2(2/ii,2/2i)), (4.1) 

iG/o *G/i iG /2 

where /o, /i, and /2 are given in (12.61) . (|2.7I) and (12.81) . respectively. 

We can obtain the maximum likelihood estimations (MLE’s) of the parameters by maxi¬ 
mizing ^ (0) in (j4.1j) with respect to the unknown parameters. This is clearly a six-dimensional 
problem. However, no explicit expressions are available for the MLE’s. We need to solve six non¬ 
linear equations simultaneously, which may not be very simple. T he maximization can be per- 
forni ed using a command like the nlminb routine in the R software (|R, Development Core Teaml . 


201J). But, it is related to initial guesses. Therefore, we present an expectation-maximization 
(EM) algorithm to find the MLE’s of parameters. 











For given n, consider independent random variables {Zi\N = n}, i = 1, 2,3 have the GLFR 


distribution with parameters nai/3 and 7 . It is well-known that 


{y^|A^ = n} = {max(Zi,Z 3 ) \N = n} , {Y 2 \N = n} = {max(Z 2 , Z 3 ) |iV = n.} 


Assumed that for the bivariate random vector there is an associated random 


vectors 


Ai = 


0 yi = Zi 

1 yi = Z 2 


and A 2 


0 y 2 = yi 

1 ^2 = ^3- 


Note that if Yi = Y 2 , then Ai = A 2 = 0. But if Yi < 12 or Yi > Y 2 , then (Ai, A 2 ) is missing. 
If (Yi, Y 2 ) € -Ii then the possible values of (Ai, A 2 ) are (1, 0) or (1,1), and If (Yi, Y 2 ) G I 2 then 
the possible values of (Ai, A 2 ) are ( 0 , 1 ) or ( 1 , 1 ) with non-zero probabilities. 

We form the conditional ‘pseudo’ log-likelihood function, conditioning on N, and then re¬ 
place N by E(N\Yi, Y 2 ). In the E-step of the EM-algorithm, we treat it as complete observation 
when they belong to Iq. If the observation belong to Ii, we form the ‘pseudo’ log-likelihood 
function by fractioning (^ 1 , 1 / 2 ) to two partially complete ‘pseudo’ observations of the form 
{yi,y 2 ,ui ( 0 )) and {yi,y 2 ,U 2 ( 0 ))) where ui ( 0 ) and U 2 ( 0 ) are the conditional probabilities 
that (Ai,A 2 ) takes values (1,0) and (1,1), respectively. Since 


we have 


PiZs < Zi < Z2\N = n) 
P{Zi < Z 3 < Z 2 \N = n) 


Ol\Oi2 


Ul ( 0 ) = 


Oi\ 


Oil -|- 03 


{oil + CK 3 ) (cKl Y Cl 2 + Ck 3 ) 
02^3 

(oi -|- Q!3) [oil Y 0.2 Y 0 ^ 3 ) 

U2 ( 0 ) = 


03 


Oil -\- 013 

Similarly, If the observation belong to I 2 , we form the ‘pseudo’ log-likelihood function of the 
from {yi,y 2 ,vi (0)) and {yi,y 2 ,V 2 (0)), where vi (0) and V 2 (0) are the conditional probabil¬ 
ities that (Ai,A 2 ) takes values (0,1) and (1,1), respectively. Therefore, 


Vl ( 0 ) 


0 X 2 

012 Y 013 ^ 


V 2 {G) 


«3 

0 X 2 + 0 x 3 


Eor brevity, we write ui ( 0 ), U 2 ( 0 ), ui ( 0 ), V 2 ( 0 ) as ui, U 2 , vi, V 2 , respectively. 

E-step: Consider bi = E{N\yii,y 2 i,&). The log-likelihood function without the additive 
constant can be written as follows: 


Veudo(0) = log( 6 ») '^bi- mlog{C{e)) + ^ log(/3 - yy*) -F ^ log(/3 - jyu) 

i=l iGIq i^IiUl2 

+ ^ log{(3 -'yy 2 i) Y {m^ui+m 2 )log{ai) Y {mi Y m 2 Vi)log{a 2 ) 
iehuh 



+{mQ + miU2 + m2V2)\og{a3) - PC^yi+ {yii + y2i)) 

i€/o iE/iU /2 

-\C^yi+^iyii + yli)) 

i€lo i€l 2 


+ai I 6 jlog(l — e 22 /i)+ 6 jlog(l — e 2 ^ 1 % 

V *€/o 


ie/iU /2 


+a 2 I ^ 6 jlog(l — e 22 ^ 0 + &ilog(l — e 2222 *^ 

I i€-fo 


iG/lU /2 


+“3 ( X] ftilog (1 - e 222 ^ + ^ 5 .l 0 g ('l 


I ie/o 


iG/i 




+ ^ 6 ilog (1 


_ p-yy2i-^yii 


2 »2i I I _ 
iG /2 / iS/i 


y^iogfi 




E 

'iG/lu /2 


_ e-/3221i-i22?i^ _ - (=“^222i-^yii 


iG/lU /2 


M-step: At this step, ^pseudo (©) is maximized with respect to the parameters. For fixed /3 
and 7 , the maximization with respect to ai, 02 and as occurs at 


di (/3,7) 
d2 (/3,7) 
as (/3,7) 


miUi+m 2 

Eie/o 7) + EiG/iu/2 (2/I0 /?, 7) ’ 

mi + m 2 Ui 

Eie/o ^**3 (2/2; 7) + EiG/iu/2 (^2*; /?, 7) ’ 

mo + mitt2 + m2U2 

T.iGio^iQiyuP,'!) + EiG/i (yii;/3>7) + Hieh biQiy2i]/3,'y) 


(4.2) 

(4.3) 

(4.4) 


where Qiy,^,^) = log(l — e 222 ^). 
solving the non-linear equation 


The maximization of f’pseudo (®) can be obtained by 


eC ( 6 ) _ - 


(4.5) 


with respect to 0 , where b = ^ Z) 2 =i bi- 

Finally, the maximization of f'pseudo (©) with respect to (3 and 7, can be obtained by max¬ 
imizing f?pseudo (di (/?, 7), d2 (/?, 7), ds (/ 3 ,7), / 3 ,7), the pseudo-profile log-likelihood function of 
13 and 7. 

The following steps can be used to compute the MLE’s of the parameters via the EM 
algorithm; 

Step 1: Take some initial value of 0, say 0^°)= ^a^°\ a^^^^ ag'^^ 7 ^°\ 0®^ . 

Step 2: compute bi = E 2 / 2 *; 0®^ 

Step 3: Compute ui, U 2 , ui, and V 2 - 



Step 4: maximize the pseudo-profile log-likelihood function ^pseudo(di(/3, 7 ), 02 (/3, 7 ), d 3 (/ 3 , 7 ), 
/ 3 , 7 ) with respect to /3 and 7 , say and 7 ^^^ respectively. 

Step 5: Compute 7 ^^^), i = 1,2,3 from (j4.2h - (l4.4p . 

Step 6: Find 9 by solving the equation ()4.5I1 . say 9^^\ 

Step 7: Replace 0^°^ by 0*^^^= 7 ^^^ , go back to step 1 and continue 

the process until convergence take place. 


5 A real example 


Che d ata set was first published in “Washington Post” and is available in 


CsorgQ and Welsh 


(|l989l f . It is obtained from the American Football League for the matches played on three 
consecutive weekends in 1986. Here, Xi represents the ‘game time’ to the first points scored by 
kicking the ball between goal posts, and represents the ‘game time’ to the first points scored 
by moving the ball into the end zone. The data are given in Table [2j 


Table 2: Scoring times (in minutes) for the matches. 



2.05 

9.05 

0.85 

3.43 

7.78 

10.57 

7.05 

2.58 

7.23 

6.85 

32.45 

8.53 

31.13 

14.58 

^2 

3.98 

9.05 

0.85 

3.43 

7.78 

14.28 

7.05 

2.58 

9.68 

34.58 

42.35 

14.57 

49.88 

20.57 


5.78 

13.80 

7.25 

4.25 

1.65 

6.42 

4.22 

15.53 

2.90 

7.02 

6.42 

8.98 

10.15 

8.87 

^2 

25.98 

49.75 

7.25 

4.25 

1.65 

15.08 

9.48 

15.53 

2.90 

7.02 

6.42 

8.98 

10.15 

8.87 


10.40 

2.98 

3.88 

0.75 

11.63 

1.38 

10.53 

12.13 

14.58 

11.82 

5.52 

19.65 

17.83 

10.85 

^2 

10.25 

2.98 

6.43 

0.75 

17.37 

1.38 

10.53 

12.13 

14.58 

11.82 

11.27 

10.7 

17.83 

38.07 


We divided all the data by 100. Then, six special cases of BGLFRPS distributions are 
considered: BGLFR, BGLFRG, BGLFRP, BGLFRB (with k = 10), BGLFRNB (with k = 2) 
and BGLFRL. Using the proposed EM algorithm, these models are fitted to the bivariate data 
set, and the MLE’s and their corresponding log-likelihood values are calculated. The results 
are given in Table [3j Eor each fitted model, the Akaike Information Criterion (AIC), the 
corrected Akaike information criterion (AICC) and the Bayesian information criterion (BIG) 
are calculated. We also obtain the Kolmogorov-Smirnov (K-S) distances with the corresponding 
p-values (in brackets) between the fitted distribution and the empirical cdf for three random 
variables Yi, Y 2 and max(yi,y 2 )- Einally, we make use the likelihood ratio test (LRT) for 
testing the BGE against other models. The statistics and the corresponding p-values are given 
in Table [3l 









Table 3: The MLE’s, log-likelihood, AIC, AICC, BIC, K-S, and LRT statistics for six sub- 


mod 


els of BGLFRPS distributions. 



Distribution 

Statistic 

BGLFR 

BGLFRG 

BGLFRP 

BGLFRB 

BGLFRNB 

BGLFRL 

di 

0.0921 

0.0605 

0.0578 

0.0597 

0.01955 

0.0675 

012 

0.5722 

0.4197 

0.3896 

0.3988 

0.1325 

0.4720 

ds 

1.1519 

0.7471 

0.7172 

0.7409 

0.2421 

0.8332 

$ 

9.6187 

12.0961 

11.4616 

11.2802 

11.6386 

12.2489 

7 

2 X 10"^ 

2 X 10"^ 

2 X 10-^ 

2 X 10"^ 

2 X 10“^ 

2 X 10“^ 

e 

— 

0.6128 

1.9930 

0.2326 

0.7186 

0.8053 

log(£) 

36.6700 

38.3625 

38.2328 

38.1661 

38.1721 

38.3582 

AIC 

-63.3400 

-64.7250 

-64.4657 

-64.3323 

-64.3443 

-64.7164 

AICC 

-61.6734 

-62.3250 

-62.0657 

-61.9323 

-61.9443 

-62.3164 

BIC 

-54.6517 

-54.2990 

-54.0396 

-53.9063 

-53.9183 

-54.2904 



0.1880 

0.1887 

0.1884 

0.1890 

0.1867 


(0.1282) 

(0.1028) 

(0.1005) 

(0.1016) 

(0.0995) 

(0.1071) 

K-S {Y 2 ) 

0.1411 

0.1469 

MhfiM 

0.1506 

0.1507 

0.1422 

(p-value) 

(0.3408) 

(0.2953) 

(0.2679) 

(0.2688) 

(0.2681) 

(0.3321) 

K-S (max(Yi,y2) ) 

0.1350 

0.1378 


0.1429 

0.1425 

0.1325 

(p-value) 

(0.3929) 

(0.3685) 

(0.3271) 

(0.3262) 

(0.3292) 

(0.4165) 

LRT 

— 

150.0651 

149.8058 

149.6724 

149.6844 

150.0565 

(p-value) 

— 

0.0000 

0.0000 

0.0000 

0.0000 

0.0000 
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